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We present the general theory of a quantum heat machine based on an A'-level system (working 
medium) whose N — 1 excited levels are degenerate, a prerequisite for steady-state interlevel co¬ 
herence. Our goal is to find out: To what extent is coherence in the working medium an asset for 
heat machines? The performance bounds of such a machine are common to (reciprocating) cycles 
that consist of consecutive strokes and continuous cycles wherein the periodically driven system is 
constantly coupled to cold and hot heat baths. Intriguingly, we find that the machine’s performance 
strongly depends on the relative orientations of the transition-dipole vectors in the system. Per¬ 
fectly aligned (parallel) transition dipoles allow for steady-state coherence effects, but also give rise 
to dark states, which hinder steady-state thermalization and thus reduce the machine’s performance. 
Similar thermodynamic properties hold for N two-level atoms conforming to the Dicke model. We 
conclude that level degeneracy, but not necessarily coherence, is a thermodynamic resource, equally 
enhancing the heat currents and the power output of the heat machine. By contrast, the efficiency 
remains unaltered by this degeneracy and adheres to the Carnot bound. 

PACS numbers: 03.65.Yz,05.70.Ln 


I. INTRODUCTION 

The extensive efforts over the years to reconcile ther¬ 
modynamics with quantum mechanics [1-5] have not yet 
fully resolved the fundamental question: What is truly 
quantum about quantum thermodynamics? Is there 
more to it than just restating traditional thermodynamic 
principles for quantized systems? Attempts to cope with 
this problem have revolved around possible quantum re¬ 
sources that may boost the thermodynamic performance 
of heat machines as compared to their classical counter¬ 
parts. A prime contender for such a resource is quan¬ 
tum coherence [6-12]. Intriguing schemes have predicted 
power and/or efficiency increase in quantum heat en¬ 
gines [7] due to sustainable (steady-state) coherence in 
non-thermal baths, as well as coherence-enhanced per¬ 
formance of photovoltaic solar heat converters [13-17] or 
coherent effects in photosynthesis [10, 18[. According to 
an interesting view [9] , coherence in thermodynamics acts 
as an information reservoir or Maxwell’s demon, i.e., as 
an extra resource that can tip the entropy balance and 
the division of input energy between heat and work in 
favor of the latter [19, 20] . 

Here we explore the thermodynamics of multilevel sys¬ 
tems with excited-state degeneracy, which is a prerequi¬ 
site for the persistence of interlevel quantum coherence 
in steady state imposed by thermal baths [6, 8, 11, 21] . 
Our goal is to address the question: To what extent is 
coherence an asset when such systems are employed as 
working media in heat machines? 

We show that the bounds on the efficiency and power 
output in the presence of coherence are general and com¬ 
mon to all cycles, i.e., reciprocating cycles that consist 
of consecutive strokes and continuous cycles. Similar 
thermodynamic performance bounds are shared by mul¬ 


tilevel systems whose excited states are degenerate and 
transition dipoles are perfectly aligned, and by multi¬ 
partite Dicke systems [22] . These general considerations 
(Sec. II) are followed (Secs. III-V) by a study of a quan¬ 
tum heat machine containing a working medium that has 
iV — 1 degenerate upper levels, whose transition dipoles 
to the ground state may not all be aligned. The system 
is constantly coupled to two spectrally distinct, hot and 
cold, thermal baths and is periodically modulated (Stark- 
shifted) by an external field. This modulation acts as a 
piston in the heat machine [23] . Our objective is to study 
the steady-state operation of such a heat machine, i.e., 
the limit cycle of its dissipative evolution, by deriving its 
heat currents, power, and efficiency in the ususal regime 
of weak system-bath coupling [1-5, 24-27] (Sec. VI). 

The present theory extends our previous study of heat 
machines whose working medium is a periodically mod¬ 
ulated two-level system (TLS) that is continuously cou¬ 
pled to two (hot and cold) baths [21, 23, 28] (see also 
Ref. [29] for a recent review). The merit of such models 
is that they are amenable to a full quantum-mechanical 
analysis by the Floquet method [23, 30, 31] . Moreover, 
their continuous-cycle operation avoids the difficulty of 
ensuring compatibility with the laws of thermodynam¬ 
ics, in contrast to machines that are operated via re¬ 
ciprocating cycles that consist of four strokes (e.g., the 
Garnot- or the Otto cycle), where, alternately, only a hot 
or a cold bath is coupled to the working medium at any 
time [32] . The difficulty therein is to properly account for 
the highly nonadiabatic energy and entropy flows induced 
by frequent on- and off-switching of interactions with al¬ 
ternate heat baths in consecutive strokes. By contrast, 
continuous-cycle, periodically-driven, heat machines can 
in a straightforward manner be made compatible with 
the first and second laws of thermodynamics regardless 


2 


of their nonadiabaticity [5, 23]. 

Steady-state coherences between upper levels are 
shown here to arise as a consequence of the thermal- 
ization of the initial state under the condition of strict 
transition-dipole alignment. However, thermalization 
may then be partially blocked, i.e., be incomplete, by the 
mere presence of dark states. For this reason we intro¬ 
duce the notion of thermalization capability of the initial 
state. We show that this capability should be maximized, 
by avoiding the initial population of dark states, in or¬ 
der to cause maximal power enhancement by the mul¬ 
tilevel heat machine, as compared to its two-level coun¬ 
terpart. We argue that it is the thermalization capabil¬ 
ity rather than steady-state coherence that underlies the 
power-enhancement mechanism: The key resource of the 
multilevel heat machine is the thermalization of transi¬ 
tions that share a common ground state, and are thereby 
correlated, whether coherently or incoherently. Similar 
principles govern heat machines based on a multiatom 
Dicke system, notwithstanding their different dynamics 
(Sec. VII). 

The heat currents and the power are analytically shown 
(Secs. VI and VII) to be both strongly boosted by the 
same enhancement factor with respect to a single two- 
level system. Consequently, the efficiency of the mul¬ 
tilevel or multipartite heat machine equals that of its 
two-level counterpart and thus adheres to the Carnot ef¬ 
ficiency bound, which is reached at zero power—the op¬ 
erating point at which the quantum heat engine is trans¬ 
formed into a quantum refrigerator [23]. The Curzon- 
Ahlborn limit [33] for the efficiency of classical heat en¬ 
gines at maximum power, however, can be surpassed for 
appropriately engineered bath spectra and carefully cho¬ 
sen temperatures [23]. 

In Sec. VIII we discuss possible realizations of the per¬ 
tinent models and their limitations. Notably, we point 
out that the iV-atom Dicke system may be the most 
straightforward experimental implementation of the ef¬ 
fects predicted in this work. In Sec. IX the conclusions 
of this work are presented. 

II. THERMODYNAMIC ROLE OF 
COHERENCES AND DARK STATES 

A. Model 

We consider a generalized V-type A-level system con¬ 
sisting of iV— 1 degenerate excited states 11),..., jiV — 1) 
and a common ground state jO). The transition frequency 
between the excited states and the ground state is de¬ 
noted by wq. This system is assumed to be coupled 
(alternately or constantly) to cold and hot heat baths 
and periodically-driven (frequency-modulated) by a pis¬ 
ton (see Fig. 1). Energy can flow from and to the baths 
via absorption and (stimulated and spontaneous) emis¬ 
sion, respectively, of thermal bath quanta. 

The system-bath interaction Hamiltonian is assumed 



Figure 1. (Color online) (a) Sketch of a heat machine. In a 
reciprocating-cycle implementation the cold and hot baths as 
well as the piston are alternately coupled to the working fluid 
during specific strokes, whereas in continuous-cycle operation 
the system interacts with both baths as well as the piston at 
all times. The sketch shows the heat machine operating as 
an engine. In the refrigeration mode all energy flows (arrows) 
have to be reversed, (b) A quantum heat machine based on 
an A-level system with {N — l)-fold excited-state degeneracy 
that constantly interacts with a cold (dashed blue) and a hot 
(solid red) heat bath. The transition frequency is periodically 
modulated to allow for power extraction (heat engine) or sup¬ 
ply (refrigerator), (c) An analogous quantum heat machine 
based on an ensemble of two-level atoms. 

to be of a generalized dipole-coupling form [25] that is 
expressed in the rotating-wave approximation as [21] 

N-l 

Hsb=Y1 Y. + , ( 1 ) 

i = l iG{c,h} 

in terms of the excitation and de-excitation Panli opera¬ 
tors 

ai-.= \j){Q\ (2a) 

(2b) 

the (possibly complex) dipole moments for the jth tran¬ 
sition dj, and the conpling operators of the cold (hot) 
heat bath Be (Bh). We adopt the notation 

aj := |di|/|di|, (3) 

where the strength of the largest transition dipole jdij 
serves as reference (i.e., aj < 1). We further define the 
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dipole alignment matrix 


:= 


d*d, 

Id, I Id, 


= cosZ(d,,dj), (4) 


where pij G [ 0 , 1 ], which encodes the directional config¬ 
uration of the — 1 dipoles. This matrix is Hermitian, 
as pij = pji and ipij = —(fji- In what follows we will for 
brevity denote dipoles d^ and dj with pij = 1 as being 
parallel even if they differ by a phase {ipij ^ 0 ). 

Prior to discussing in detail a possible setup of such a 
heat machine, we infer from the system-bath interaction 
Hamiltonian (1) some general results concerning the ther¬ 
modynamic role of coherences in the working medium. 


If, however, the initial state of the multilevel system 
is non-dark, then the rate of quanta exchange of the ef¬ 
fective TLS (formed by Itjjh) and |0)) with the baths, 
7 b, is enhanced by a factor of — 1 compared to the 
spontaneous-emission rate 71 of a TLS consisting of the 
states |0) and |1). This can be deduced from the (dissipa¬ 
tive) master-equation description for the reduced density 
matrix of the system, based on the interaction Hamilto¬ 
nian (7): The spontaneous-emission (decay) rate scales 
with the square of the transition-dipole moment [24] and 
therefore 7 b = (A^ — l) 7 i with 71 oc |dip (Sec. VB). 


C. Power and efficiency bounds 


B. Collective-states basis 

Owing to the degeneracy of the excited states of the 
system Hamiltonian 

N-l 

Hs = cr^+crf, (5) 

i=i 

one is free to choose any (rotated) basis within this ex¬ 
cited manifold, as the rotated states will still be energy 
eigenstates of the Hamiltonian (5). Assuming that all 
dipoles are parallel and (for simplicity) of equal strength 
and real, we then transform the system-bath Hamilto¬ 
nian ( 1 ) to the basis consisting of the ground state | 0 ), 
the collective “bright” state [34] 

= (6) 

and N — 2 rotated excited states that are orthogonal to 
the latter. The system-bath interaction Hamiltonian (1) 
then adopts the form 


Let us consider the implications of the collective en¬ 
hancement for the heat currents, i.e., the rate of en¬ 
ergy exchange between the system and the baths, and 
the power of the heat machine. The heat currents 
must be proportional to the rate of quanta absorption 
and emission of quanta from and to the baths, respec¬ 
tively [3, 4, 36]. This proportionality should be indepen¬ 
dent of the actual implementation of the heat machine, 
i.e., whether both baths are continuously coupled to the 
system or only one at a time in a reciprocating cycle. 

Due to the enhancement of the decay rate for parallel 
dipoles in a bright state, the “cold” (Jc) and the “hot” 
(Jh) heat currents Ji = Qi, Qi being the heat exchange 
with the fth bath, will be equally enhanced compared to 
their two-level counterparts, as 

Je = (Af-l)J™ ( 8 a) 

Jh = (fV-l)Jh^L®. ( 8 b) 

The first law of thermodynamics (energy conserva¬ 
tion) [3, 0 , 23] then implies that the power is enhanced 

by the same factor, 

IT = -(J^-HJh) = (Af-l)IT™. (9) 


^parallel ^ ^ (IV'b) (0| 0 di • B, + H.C.) . 

iG{c,h} 

( 7 ) 

Formally, this Hamiltonian describes a single two-level 
system formed by | 0 ) and |' 0 b)j whose interaction with 
the baths has a dipole moment enhaneed by a factor of 
y/N — 1, which is responsible for superradiance [22, 35] . 
The remaining N — 2 states, which are orthogonal to |0) 
and I'lph), are not accessible by this Hamiltonian: They 
are dark states with respect to the (dipolar) system-bath 
interaction. Consequently, if the A^-level system with 
parallel dipoles is initially prepared in one of these dark 
states, it does not exchange heat with the baths, and its 
state remains invariant under the action of the system- 
bath interaction Hamiltonian. We therefore anticipate 
the result (elaborated further on), that the steady-state 
solution for parallel dipoles strongly depends on the over¬ 
lap of the initial state p( 0 ) with these dark states, i.e., 
on the initial value (nd)p(g) of the dark-state projector. 


Equations ( 8 ) and (9) imply that although an enhance¬ 
ment of the output power is expected, the efficiency of 
the multilevel heat engine, which is defined by the ratio 
of the extracted power to the invested “hot” current [32] , 


V ■= 


A 


|t^TLS| 




^TLS ^ 


( 10 ) 


remains the same as for a TLS-based heat engine where 
coherence is absent. Notably, the Carnot bound [32] is 
adhered to, by virtue of Eq. (9), Eq. (10) and the second- 
law condition [37] 


A 

A 



( 11 ) 


The assumptions in the foregoing discussion were the 
degeneracy of the excited states as well as the dipolar- 
coupling form of L7 sb between the system and both ther¬ 
mal baths. Any increase of the efficiency compared to 
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Eq. (10) would require different enhancement of Jc and 
Jh, such that the enhancement Ch of the “hot” current 
J\i = exceeds its “cold” current counterpart Cc, 

i.e., Ch > Cc- This may only be possible for different 
coupling Hamiltonians between the system and the two 
thermal baths, and HgB, unlike Eq. (1). However, 
the Carnot bound, which is a corollary of the first and 
second laws [32], is upheld regardless of these assump¬ 
tions since even if such enhancement can be achieved, 
the second law restricts the ratio of the current enhance¬ 
ment factors Cc and Ch to satisfy [cf. Eq. (11)] 

rTLS tTLS i 

"'c , ^ n 



D. Role of coherences 

If all dipoles are parallel, coherences in the working 
fluid (medium) will build up in the bare (energy) basis 
according to the Hamiltonian (7) once the system is in 
contact with heat baths at finite temperature, since the 
bright state IV'b) [see Eq. (6)] is a coherent superposition 
of the bare excited states. Initial coherences may either 
suppress or enhance the machine performance: They can 
bring the heat machine to a complete standstill if they 
correspond to the set of dark states liv'd)} since the 
Hamiltonian (7) then yields iJgB |V'd) ~ 0, so that the 
system is then decoupled from the baths. Conversely, 
the Hamiltonian (7) will enhance both the heat currents 
and the power by a factor of N— 1 due to the superradiant 
(collective-decay) effect if either |V’b) or the ground state 
|0) is initially prepared. Yet, even if the system has ini¬ 
tially no coherence, it will evolve via the Hamiltonian (7) 
to a thermal mixture of |0) and |V'b) that has persistent 
coherences. Thus, steady-state coherences (embodied by 
the bright state) are formed in the bare energy basis as a 
result of thermalization, but the entire system does not 
thermalize, i.e., an iV-level Gibbs state does not arise be¬ 
cause of the availability of dark states (even if they are 
not populated) (see Sec. VB). 

The foregoing discussion of the thermodynamic role of 
coherences in the working medium holds for the special 
case of all transition-dipole vectors in the Hamiltonian (1) 
being parallel. In what follows, we, however, wish to ex¬ 
plore the possibility of power enhancement for misaligned 
(non-parallel) dipoles, or any combination of aligned and 
misaligned dipoles. Since steady-state coherences corre¬ 
spond to (destructive or constructive) interference of the 
dipoles, they can only arise if at least some of the dipoles 
are aligned. Only then may dark states exist. With¬ 
out the presence of such dark states, a thermal (Gibbs) 
steady state of the working fluid is expected. 

As we show in Sec. VH, these arguments also hold for 


an ensemble of N two-level atoms in a suitable geometry, 
which realizes the Dicke model [22, 35] (Sec. VHIG). In 
that case the effective system-bath interaction Hamilto¬ 
nian cannot be written in the form of Eq. (7). Instead, 
the ensemble can then be mapped onto a collective spin- 
N/2 system and a multitude of dark states. The heat 
currents will no longer adhere to a simple form as in 
Eqs. (8), but the efficiency of the Dicke heat machine 
will still be the same as of the TLS-based machine. 


III. FLOQUET EXPANSION OF THE MASTER 
EQUATION IN CONTINUOUS CYCLES 

In the remainder of this article these issues are inves¬ 
tigated for a continuous-cycle heat machine, wherein the 
{N — l)-fold degenerate system is constantly coupled to 
the two thermal baths and to a periodically modulating 
“piston” that allows for work extraction (in the case of 
an engine) or supply (in the case of a refrigerator or heat 
pump), respectively (see Fig. 1). 

In our model the “piston” effects are described by 
a synchronous periodic modulation of all the upper 
states [21, 23, 28], 

N-l 

Hs{t) = h[uJo+OJit)] ^ cr|_cri, (13) 

i=i 

where uj{t + ^) = w(t), Q being the modulation rate. 
Such a transition-energy modulation may for example be 
induced by a varying magnetic field (via the Zeeman ef¬ 
fect) or by an alternating electric field (via the Stark 
effect). 

The master equation [24] for the reduced density op¬ 
erator p in the interaction picture is 

P = (14) 

where the Liouvillian superoperator £ is of the Lindblad- 
Gorini-Kossakowski-Sudarshan (LGKS) form [26, 27]. 
As shown in [23, 30] this Liouvillian can be decomposed 
into “sub-bath” Liouvillians £| describing the interaction 
of the system with the ith bath {i G {c,h}) evaluated at 
the gth harmonic sideband (q G Z) of the unperturbed 
transition frequency wq, induced by the modulation rate 
Q. These “sub-bath” (sideband) contributions are a con¬ 
sequence of the Floquet theorem for the solution of lin¬ 
ear differential equations with periodic coefficients. The 
master equation (14) then adopts the additive form 

P = Y1 H 

g£Z i—{c,h} 

with the “sub-bath” Liouvillians (generalizing Ref. [21]) 
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N-l 


C‘}p=lpiq)G.{uJo + qn)Y, 


i=i 


a'^j'D (cri,cr+) + pjye*‘^«'ajajvX> 


N-l 


+ 7^P{<l)Gi{—uJo — q^) 'Y 


i=i 


aji? (o'+^o--) + Y "'^"'ajafV (<,cri) 
iVi 


. (15b) 


Here the dissipator is defined as 'D{a,b) := 2apb — 
hap — pba for arbitrary system operators a,b. The terms 
P (^cr^-,cr^+^ and V ^cr;^, in Eq. (15b) describe (spon¬ 
taneous and stimulated) emission into and absorption 
from the bath, respectively, involving a single transition 
dipole dj and weighted by [cf. Eq. (3)]. By virtue of 
sharing a common ground state |0), these contributions 
to the Liouvillian amount to an indirect population trans¬ 
fer between the different excited states via this common 
ground state. 

By contrast, the cross-terms and 

P in Eq. (15b), weighted by ajaj', describe 

correlated absorption and emission involving two differ¬ 
ent transitions, i.e., quanta exchange between excited 
states \j) and \j') via the common ground state |0), 
\j) —?> |0) —)■ \j'). These bath-mediated interactions 
between the excited states—which result in dynamical 
coherences—are largest when the two dipole moments 
involved, and d^/, are parallel up to a phase fac¬ 
tor (pjj/ = 1). For orthogonal dipole orientations 
{pjji = 0 Vj 7 ^ j') these cross-correlations vanish and no 
coherences build up. In three-dimensional space, this 
completely orthogonal configuration can only be realized 
in a three- or four-level system. 

The rates of all decay and absorption processes are 
determined by the respective prefactors of the dissipa- 
tors in the Liouvillian (15b). Here P{q) are the Floquet 
coefficients [30] determining the weight of the gth har¬ 
monic sideband (with normalization Pi'l) ~ 1) 

Gi{uj) is the response spectrum of the ith bath evaluated 
at frequency w. These spectra fulfill the Kubo-Martin- 
Schwinger (KMS) detailed-balance condition [25] 

G,{-uj) = (16) 

For a bosonic bath, the coupling strengths explicitly read 

Gi(w) = 7i(w) (ni(a;)-I-1), (17) 

with ni{uj) := — 1)“^ denoting the number of ther¬ 

mal quanta at inverse temperature fdi = l/k^Ti and 7 i(u;) 
being the frequency-dependent transition rate induced by 
the ith bath. 


IV. THREE-LEVEL SYSTEM 

In Sec. HI we have introduced the full model, invok¬ 
ing couplings to two heat baths and periodic modulation 
of V — 1 excited states. These ingredients are neces¬ 
sary for the operation of the system as a heat machine. 
However, all the “sub-bath” Liouvillians Cf [Eq. (15b)[ in 
the master equation (15a) contain the same dissipators 
P and only differ by their respective prefactors. Owing 
to this additive structure, in this section we first solve 
the master equation for a three-level system interacting 
with a single bath at inverse temperature /? = I/ZcbT 
and a static (unmodulated) transition frequency {q = 0, 
P(0) = 1), and analyze its steady-state solution. This 
solution will serve to obtain the “global” solution for 
two baths involving a modulated transition frequency 
by simple changes of the coefficients (i.e., the transi¬ 
tion rates). The two-bath “global” solution allowing for 
the modulation-induced harmonic sidebands, which is re¬ 
quired for the computation of the heat currents, will be 
presented in Sec. VI. 


A. Steady-state solution for degenerate excited 
states 

In order to obtain some insight into the dynamics, we 
transform the master equation (15a) into a set of Ehren- 
fest equations of motion for operator expectation values. 
These equations can subsequently be cast into an inho¬ 
mogeneous linear ordinary differential equation (ODE) 
system 

X = -I-b (18) 

for the density matrix elements. In the case N = 3 this 
vector is 

X := {p 2 i,Pi 2 , Poo, P 22 )^■ (19) 

The corresponding coefficient matrix A and the inhomo¬ 
geneity b are presented in Appendix A. 

This linear system does not describe the entire density 
matrix, as it does not contain the coherences poi ^^nd po 2 
between the ground and the excited states. The reason 
is that those matrix elements are completely decoupled 
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from X and obey the independent (homogeneons) differ¬ 
ential equation 

y = By (20) 

with 

y'■= iPl0,P01,P20,P02)'^ (21) 

and the coefficient matrix B given in Appendix A. All 
eigenvalues of the matrix B are negative. Consequently, 
the coherences between the ground and excited states 
are damped out regardless of their initial values and the 
steady-state solution of Eq. (20) has no coherences, i.e., 
y’*’' = 0. 


p=l p=l 



1 


’* 












r\ : 





01234560123456 


Let US now revisit the ODE (18) for x. The uniqueness 
of its steady-state solution depends on the determinant 
of the coefficient matrix A, 


det(A) = 


1 4 


2G(a;o) 


(l -F 2e~P^°) x 

+ ( 22 ) 


For aligned dipoles (p := pi 2 = 1) we find det(A) = 
0, corresponding to a singularity of the coefficient ma¬ 
trix (Al). In this regime, multiple steady-state solutions 
of Eq. (18) may exist (depending on the initial condi¬ 
tions), each satisfying the linear system of equations 


Ax"" = -b. 


(23) 


The general steady-state solution of the linear 
ODE (18) (Eq. (A5) in Appendix A) is compatible with 
the one found in Refs. [6, 38]. While the steady state is 
diagonal and unique for non-aligned dipoles (p 1), it 
depends on the initial conditions (via an integral of mo¬ 
tion) and yields persistent coherences = (pli)* the 
aligned case (p = 1). The integral of motion agrees with 
a previous finding (involving a zero-temperature bath 
with [8] and without [38] incoherent external driving, 
respectively). Such initial-condition dependent steady- 
state solutions have also been found for degenerate A- 
systems [39]. 

The time evolution of the system towards the steady 
state (A5) is illustrated in Figs. 2 and 3, where we have 
numerically integrated the master equation (15a) and 
compared its steady state to the analytic steady-state 
solution (A5). The time evolution of the ODE (18) 
gives the same result. The dependence of the steady- 
state solution on the initial conditions for aligned dipoles 
(p = 1) can clearly be seen. Coherences between the 
excited states persist only in this regime and are just a 
transient effect for any misalignment. Notably, as seen 
from Eqs. (A5), the steady states for orthogonal dipoles 
(without transient coherences) and any other misalign¬ 
ment coincide. Consequently, dynamical correlations be¬ 
tween the upper states (induced by the cross-terms in the 
Liouvillian (15b)] are important at finite times but not 
in the long-time limit. In particular, this implies that 
coherences cannot play any thermodynamic role in mis¬ 
aligned dipole transitions in the steady-state operation 
of the heat machine. 


yi -yt 

Figure 2. (Color online) Time evolution of the popula¬ 
tions (solid blue: p 22 , dashed red: poo) for aligned (p = 1 ) 
and misaligned (p = 0.7) dipoles and different initial con¬ 
ditions p( 0 ) = | 0 ){ 0 | (left) and p( 0 ) = \'tp)('tp\ with l"!^) = 
(|0) + |l)-|2))/y3 (right) obtained by numerically integrat¬ 
ing the master equation (15a). The horizontal solid black lines 
are the analytic steady-state solution (A5) of the ODE (18). 
Parameters: = |, a = 1 and G{ujo) = 'y{n -I- 1) = 27 . 


p = l 


P = 


1 



p = 0.7 p = 0.7 



01234560123456 
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Figure 3. (Color online) Same as Fig. 2 for the coherence P21. 
In the aligned case (p = 1) the coherence persists in steady 
state, whereas it is merely a transient effect in the misaligned 
case (p = 0.7). 


B. Non-degeneracy effects 


The singularity of the coefficient matrix A for aligned 
dipoles only arises for perfectly degenerate levels. If a 
detuning A between the two excited states is introduced 
into the coefficient matrix (Al), the determinant becomes 
proportional to [(1 -|- Q;^)^(l — p^) -f A^/G^(a;o)] in the 
Schrodinger picture ]as a generalization of Eq. (22)]. 
Amongst the solutions of det(A) = 0, only the case 
p = 1 and A = 0 lies within the domain of existence 
of p, which is [0,1]. In the interaction picture (in which 
the master equation (15a) is presented) the coefficients 
then become time dependent through the phase factors 
pe^*”^ I-! and the steady state is no longer 
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determined by the determinant of the coefficient ma¬ 
trix. In the secular approximation, these time-dependent 
phase factors, and hence all cross-terms in the Liouvil- 
lian (15b), vanish upon coarse-graining the time evolu¬ 
tion [25]. The steady-state solution is then given by the 
canonical (Gibbs) distribution (diagonal in the energy 
basis) plf = e~^^'pQ^Q for 1 < i < N—1. Hence, we must 
assume degenerate levels in order to study possible im¬ 
plications of quantum coherences on the thermodynamic 
properties of a steady-state quantum heat machine. 


B. Parallel dipoles 

Let us now consider the opposite case of all dipoles 
being parallel (which has been previously discussed in 
Sec. II). The steady-state solution in this case may not 
be unique due to the singularity of the coefficient matrix, 
which is linked to the existence of dark and bright states. 

When all dipoles are aligned with the reference dipole 
di, i.e., dj = the interaction Hamiltonian (1) 

between the system and the bath can be recast into the 
generic form [cf. also Eq. (7)] 


V. STEADY-STATE SOLUTION: DEPENDENCE 
OF THERMALIZATION ON THE DIPOLE 
ORIENTATIONS 

In the previous section we have shown that in the sim¬ 
plest case {N = 3), aligned dipoles are dynamically dis¬ 
tinct from all other possible geometric configurations (cf. 
also Ref. [21]). In this section we show that this distinc¬ 
tion arises from the existence of a dark state (anticipated 
in Sec. H), which hinders full thermalization under cer¬ 
tain initial conditions if the dipoles are all aligned. We 
expect similar behavior for the general case of A — 1 
dipole transitions. In the following we investigate the 
steady-state solution for arbitrary dipole orientations, 
taking into account that as N grows, so does the dipole 
configuration space. To this end we first consider the ex¬ 
treme, simple cases of (i) no dipoles being pairwise paral¬ 
lel and (ii) all dipoles being parallel. The general solution 
for arbitrary dipole orientations can then be derived in a 
transparent way from these extreme cases. 


A. No pairwise parallel dipoles 


If no dipoles are parallel, the master equation (15a) 
(for a single bath at inverse temperature (3 and without 
modulation, q = 0) possesses a unique thermalized, diag¬ 
onal, steady state 


= jvTTTPSiS ^■«[! 

ss ^ 1 

^00 j _ iy-phu,o ■ 


iV- 1] 
(24a) 
(24b) 


The straightforward proof is given in Appendix B. Here 
too, the steady state is independent of the actual orienta¬ 
tion of the dipoles (as long as they are not aligned)—but 
the time evolution leading to this steady state depends 
on their orientation. 

Owing to their degeneracy, any superposition of ex¬ 
cited states is again an energy eigenstate. This means 
that the density matrix (24) is not only diagonal in the 
“bare” basis {|0), |1),..., | A — 1)}, but also in any ba¬ 
sis {|0), |1'),..., |(iV — 1)')} formed by the ground state 
and the rotated excited states. 


N-l 


^^Sb| 


Pij=i 




((7+ O di • B -f cr_ (g) d* • B^) . 


(25) 

This operator is reported in the basis spanned by the 
following states: (i) the ground state |0), (ii) the bright 
state [see also Ref. [34] and Eq. (6)[ 


IV'b) := , \j ), (26) 

VE,=1 aj 

and (iii) N — 2 dark states liv’d)}’ obtained by the 
Gram-Schmidt orthogonalization procedure. 

The Hamiltonian (25) formally looks like the coupling 
Hamiltonian of a single two-level system interacting with 
the environment. However, this operator acts on an N- 
dimensional Hilbert space, even if only two levels (the 
ground and the bright states) are explicitly involved, via 
the Pauli excitation and de-excitation operators [cf. also 
Eq. (7)[ 


^+:=|V’b)(0| (27a) 

tr- := |0)(V’b|. (27b) 

The associated transition-dipole vector is scaled by a fac¬ 
tor of which for equal transition strengths 

(aj = IVj) amounts to \/N — 1. The diagonal steady- 
state solution of the master equation following from 
the interaction Hamiltonian (25) has the form (cf. Ap¬ 
pendix G) 


Phh — 


Poo — 


1 

X _l_ gpnujo 

1 

1 -|- 


1 - (Hd) 


p(o) 




1 - (Hd) 


p(o) 


PTk = (V>/c|p(0)|V’fe) for fc = 1,..., A - 2. 


(28a) 

(28b) 

(28c) 


Here we have defined the projector 


Hd := 


N-2 

E 


V’d 


= l-|V'b)(V’b|-|0)(0| (29) 


onto the dark-state subspace, the dark-state populations 
pI’)., and the spectral decomposition of the part of the 

















initial state lying within the dark-state subspace, 


N-2 N-2 

pd{o)= ^d)(v'd = ■ 

i,j—l k—1 


(30) 

This density matrix is not normalized to unity but to 
Trpd(O) = (nd)^(o) = 1 - (V’b|/o(0)|^b) - (0|p(0)|0). The 
states \(j)k) associated to the dark-state populations 
are thus the eigenvectors of pd(0). If Pd(0) = 0, any basis 
of the dark subspace may be chosen instead. 

The dark-state populations in the steady-state so¬ 
lution (28) are the same as their initial values, since these 
states do not interact with the environment. The bright- 
and ground-state components (which do interact with 
the environment) eventually thermalize, their Boltzmann 
factor being determined by the bath temperature. 

However, the steady-state thermalization is only partial 
if the initial state has an overlap with the dark-state sub¬ 
space. We therefore refer to the factor 1— (nd)p(Q) as the 
thermalization capability of the initial state p(0). 

The foregoing discussion has shown that whilst the in¬ 
terpretation of the initial-condition-dependent terms in 
the non-diagonal density matrix (A5) may be obscure, 
their physical meaning is clearly revealed in the diagonal 
steady-state form (28), which has a distinct non-Gibbs 
character. 



Figure 4. (Color online) A possible configuration exhibiting 
a domain structure of parallel dipole transitions. Each group 
(domain) of aligned transitions (upper arrows) is equivalent 
to a single (enhanced) dipole transition (lower arrows). 


\Pi) \Pl) K) |iVp + i> |iv^ + 2 > ijv-i> 



Figure 5. (Color online) Effective A^eff-level scheme corre¬ 
sponding to the A^-level system in Fig. 4. The transitions 
associated to a bright state (i.e., to a domain of parallel 
dipoles) are enhanced compared to their bare-state counter¬ 
parts. Within this effective system no transition dipoles are 
parallel to each other. 


C. Arbitrary dipole orientations 

Based on the two limiting cases considered above, we 
expect the dynamics to be different for aligned (parallel) 
and misaligned transition dipoles. 

We therefore group all aligned transition dipoles (up 
to a phase factor) into “domains”: The first ni aligned 
dipoles are grouped in domain 1, the next n 2 parallel 
ones are grouped in domain 2 and so forth. Altogether, 
there are p “domains”, which contain 

p 

(31) 

i=i 

transition dipoles. The remaining iV — 1 — Np transition 
dipoles are non-parallel to any other dipole (see Fig. 4). 

As detailed in Appendix D, the Wlevel system with ar¬ 
bitrary dipoles can be mapped onto an Wff-level system 
with non-parallel dipoles, where 

N,s-.= p + N- Np. (32) 

Upon introducing the new coupling constants 

for m < p 

aArp+(m_p) for p -I- 1 < m < Wff - 1 

(33) 


the interaction Hamiltonian adopts the form 

-/Veff — 1 

Hsb= Y. ^mldiKa^Ge^-B + d^GeJ^-Bt). 

m—1 

(34) 

This is simply the Hamiltonian (1) for levels with 
modified coupling strengths and no pairwise parallel 
dipoles. The new Pauli operators If™ are [cf. Eqs. (2) 
and (D2)] 

CT™ = IV'™) (0| m < p 

P+l<m<iVeff-l' 

= |Ap-b(m-p))(0| 

(35) 

The mapping from N levels with arbitrary transition 
dipoles to fVeff non-aligned dipoles (cf. Figs. 4 and 5) 
now assigns a weight Aj to each transition, according 
to the number of thermalization pathways forming each 
domain of parallel dipoles. Clearly, the completely mis¬ 
aligned (Wff = N) and aligned (Wff = 2) cases discussed 
in Secs. VA and VB above are included in this general 
scheme. 

The Liouvillian associated with the Hamiltonian (34) 
is thus the same as for the non-aligned system obtained 
upon replacing N i— Nes and aj i— Aj. The steady-state 
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solution of the master equation then reads 


in terms of the harmonic sideband contributions 

Poo = 

1 


1 (nd)p(o) 

(36a) 

■//:=-|Tr[(£fp-)ln(p«)]. (40) 

Pi 

1 + (Aeff - l)e 

—(ShiijQ 

/o’*® - 


i = 1, 

• ■ , Neff — 1 

(36b) 

Here p®® is the global steady-state solution of the master 

rii 


i = NeS -k 1,..., A 

equation (15a) satisfying 


Here again the {\(t>i)} are the eigenvectors of pd(0), i.e., 
the part of the initial density matrix lying within any of 
the p dark-state subspaces. The projector Hd onto the 
union of all dark-space subspaces has now the form 


P rik-l 

nd = ^ |V'd'')(V’d^ 

j=l k=l 


(37) 


N-1 


= li)0’l - |0)(0| • (38) 


i=i 


i—-^p +1 


The proof of Eqs. (36) is analogous to the respective 
proofs of the steady-state solutions (24) and (28) for non- 
aligned and aligned dipoles. 

Whilst the steady-state density matrix (36) is diago¬ 
nal in the new basis associated with the Hamiltonian (34) 
(and hence does not contain any coherences), off-diagonal 
elements do appear in p®®, the representation of the den¬ 
sity matrix p®® in the energy (bare-state) basis [cf. also 
Eq. (AS)]. 

It is important to note that the coherences with the 
largest modulus, arise for a dark ini¬ 

tial state, since the exchange between the excited states 
[described by Eqs. (15)] implies that the ground state 
is at least partly (thermally) populated in steady state 
[cf. Eq. (36a) [. This ground-state population inevitably 
reduces the excited-state population and hence the coher¬ 
ences between excited states. Equivalently, the steady- 
state density matrix can only describe a completely 
population-inverted state iff the initial state is dark. Yet, 
dark states are a severe hindrance to the operation of a 
quantum heat machine. Hence, the magnitude of steady- 
state coherences cannot be a criterion for possible power 
enhancement in a multilevel quantum heat engine, as 
shown below in detail. 


= 0, (41) 

and pf is the local solution for the system interacting 
with a single (ith) bath evaluated at the gth harmonic 
sideband, i.e., the steady-state solution of the Liouvil- 
lian (15b) 


= 0. (42) 

The heat currents (40) are based on a global two-bath 
solution, as required to avoid violation of the second 
law [40[. 

The steady-state solution of the entire master equa¬ 
tion (15a) can be obtained from Eq. (36) upon adapting 
the transition rates to a two-bath situation. Owing to the 
KMS relation (16), this amounts to replacing the real in¬ 
verse temperature /3 by an effective inverse temperature 
/3eff, which is defined by the Boltzmann factor 

Sketch} ^(g)G.(-^o-g») 

This effective temperature can be controlled (modified) 
by engineering the modulation type [and thus the Floquet 
coefficients P{q)] or the cold- and hot-bath response spec¬ 
tra [Gc(wo + 9 O) and Gh(wo + 9 ^)] at the corresponding 
harmonic sidebands g G Z, respectively. Likewise, the 
local steady-state solution pj is obtained from Eq. (36) 
upon replacing 


/3 !->• /3i (44a) 

ojf) I —y ojq gO. (44b) 

Inserting p®® and p| into expression (40), we find the 
heat currents, whose explicit form is presented in Ap¬ 
pendix E. Their comparison to the previously found re¬ 
sult for a TLS [23[ yields the heat current ratio 


VI. HEAT CURRENTS AND POWER 
A. Heat currents 

The general steady-state solution (36) of the master 
equation (15a) is the key ingredient for computing the 
heat flows (currents) between the system and the cold 
and hot baths. Upon invoking the dynamical version of 
the second law of thermodynamics and the von Neumann 
entropy, the heat currents associated to the Ah bath are 
found to be [5, 23, 30[ 

Jz = ^ J/ (39) 

qGlj 



'N-1 

E< 


1 (nd)p(o) 


2^ _|_ g-/3ott?ii^o 
1 -k (Veff - ■ 


(45) 

Here the first factor on the right-hand side is the enhance¬ 
ment that stems from the iV — 1 available thermalization 
pathways. It amounts to A — 1 if all transitions are of the 
same strength {aj = IVj). The second factor measures 
the thermalization capability of the initial state (i.e., the 
overlap of the initial state with the non-dark states, viz. 
the part of the initial state amenable to thermalization). 
This factor is only relevant for the case of at least two 
dipoles being parallel (for other configurations there are 
no dark states). The last factor strongly depends on the 












10 


modulation type and the bath properties via the inverse 
effective temperature PeS defined in Eq. (43), but also 
on the effective number of misaligned dipoles [cf. 
Eq. (32)]. This last factor becomes unity for Nes = 2. 
Such an effective two-level system is obtained for an N- 
level system with all transition dipoles being aligned [cf. 
Eqs. (28)]. 

Remarkably, for any number of misaligned dipoles the 
product of the last two factors in Eq. (45) is determined 
by the ratio of the respective ground-state populations, 
yielding 

j /n- 1 \ ss 

j^= ( (46) 

It follows from Eq. (46) that the largest possible en¬ 
hancement factor is bounded by 

N-l 

(47) 

i=i 

the equality sign corresponding to all transitions having 
the same strength. This maximal enhancement corre¬ 
sponds to the combined heat currents from — 1 inde¬ 
pendent two-level heat machines. 

The most advantageous initial condition is a state or¬ 
thogonal to the dark subspace, e.g., p(0) = |0)(0|. In that 
case, as is apparent from Fig. 6, at low effective temper¬ 
atures all dipole orientations yield the same (maximum) 
enhancement by a factor of iV — 1 (assuming all dipoles 
are equally strong). This can be understood from the 
population ratio in Eq. (46): At low temperatures the 
ground-state population is close to one, independent of 
the dipole configuration. 

In the opposite case of high effective temperatures, 
configurations with more parallel transitions are more 
favorable. This can again be explained based on the 
ground-state population ratio in Eq. (46). Assuming 
full thermalization, the ground state is the more pop¬ 
ulated the fewer levels are available. Hence, for N levels 
the smallest N^s are the most beneficial at such tem¬ 
peratures. For the 10-level system shown in Fig. 6, the 
ground-state population in the high-temperature limit 
(/3eff 0) is [according to Eqs. (36)] = 1/A^eff, 

compared with pS"® = 1/2, so that Eq. (46) yields 
limp^,,^o J,/JTLS = 9 X 2/lVeff. 


B. Power and efficiency 


The first law (energy conservation) W = — (Jh + Jc) 
relates the heat currents to the power. Whilst the explicit 
form of the power is again reported in Appendix E, we 
here compare it to its TLS counterpart. 


W 

fETLS 



Poo 


J^ 




<N-l. 
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Figure 6. (Color online) Heat current ratio Ji/ [Eq. (45)] 
for a ten-level system assuming no initial overlap with any 
dark states for (1) Aeg = A = 10 (no parallel dipoles), (ii) 
Aeff = 6 (some dipoles parallel) and (iii) Aefi = 2 (all nine 
dipoles parallel). Alignment is only beneficial at high effective 
temperatures (small ffshLOo). For all dipoles being parallel, 
amplification is independent on the bath temperatures. Due 
to Eq. (48) the same holds for the power ratio. 


This is the same ratio as in Eq. (46) for the heat cur¬ 
rents and therefore Fig. 6 also holds for the power ratio. 
We have thus come to a central conclusion: An A-level 
quantum heat machine can never outperform N — 1 in¬ 
dependent two-level heat machines. 

In contrast to the power, the efficiency 

-W 

r? = — (49) 

(when the heat machine acts as an engine) and the coef¬ 
ficient of performance 

COP = i (50) 

W 

(when the heat machine acts as a refrigerator), respec¬ 
tively, are not affected by the presence of degenerate up¬ 
per states, regardless of the dipole orientation. Notably, 
the Carnot bound is not violated, based on the results 
previously obtained for the TLS heat machine [23], as 
discussed in Sec. IIC. 


C. Power-dependence on modulation rate 

The result (E2) for the power output of the 
periodically-driven continuous multilevel quantum heat 
machine only differs by the extra factor JVeg — 1 in the 
denominator from its counterpart for a two-level machine 
discussed in Ref. [23]. It was shown there that depending 
on the modulation frequency the heat machine can be 
operated either as an engine (W < 0 for ft < flcrit) or 
as a refrigerator/heat pump {W > 0 for > flcrit), with 
flcrit being some critical frequency. 

This behavior can be understood in terms of the fre¬ 
quency dependence of the bath populations nffuj) and 
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nh(w). At the TLS resonance freqnency, they fulfill 
nh(wo) > nc(uJo)- The periodic driving introduces har¬ 
monic sidebands, that shift the frequencies at which the 
TLS couples to the baths [cf. also the Liouvillian (15b)]. 
When only two sidebands wq ± contribute and the 
cold and hot bath spectra do not overlap, the “nat¬ 
ural” direction of the heat flows is reversed provided 
that nh(oJo -|- fl) < nc(oJo — fl) [31]. This explains 
the occurrence of a “critical” frequency flcrit) defined by 
nh(wo -l- flcrit) = nc{uJo ~ ^^crit)j ^t which the heat flows 
(and hence also the power) vanish. This is the point at 
which Carnot efficiency is reached [23]. Similar scenar¬ 
ios occur for example for ultracold atoms in cavities ]41], 
where, depending on the frequency of the cavity-pump 
laser, the atoms are either cooled or heated. 


D. Power-dependence on initial conditions for 
aligned dipoles 

As we have seen, the power enhancement of the iV-level 
heat machine compared to a single two-level heat ma¬ 
chine strongly depends on the initial conditions if some 
of the transition-dipole vectors are aligned ]Eq. (E2)]. We 
here illustrate this dependence for a three-level {N = 3) 
system, which allows for a simple graphical representa¬ 
tion. We shall now specify its advantageous initial condi¬ 
tions. Once the ground-state population /Ooo(O) is fixed, 
the advantageous states, amenable to full thermalization, 
are those having the largest possible modulus of the co¬ 
herence and the appropriate phase. These “optimal” co¬ 
herences and populations must satisfy (cf. Sec. IV A) 

P2i(0) = ^^[l-Poo(0)] (51a) 

Pii(O) = 1 2 [^~Poo(0)] (51b) 

P22(0) = ^ [l-poo(O)]- (51c) 

1 + 

These conditions may be satisfied by a pure state p(0) = 
|^/>(0))(V^(0)|, where 

IV'(O)) = Vpoo(O) |0) + x/ 1 - poo(0)e*‘^ jV'b) , (52) 


power enhancement factor 



2 

1.8 

1.6 

1.4 

1.2 

1 

Reduction 


Figure 7. (Color online) Power enhancement factor (55) for a 
three-level system with respect to a TLS as a function of the 
ground-state population and the coherence for parallel equal 
dipole moments. The optimal amplification (by a factor of 
two) is realized on the upper edge of the colored area, i.e., for 
the maximally allowed coherence (and correct phase) once poo 
is fixed ]Eq. (51a)]. The lower left shaded triangle corresponds 
to a reduction of the power relative to a TLS, up to complete 
suppression in its lower left corner (the point corresponding 
to the dark state, which has the maximally allowed value for 
the modulus of the coherence but the opposite phase). For 
anti-parallel {ip = tt) dipoles the diagram would be vertically 
flipped around the Rep 2 i = 0 axis. 


reduced power output. The power-enhancement factor 
for two equal parallel dipoles (a = 1 and tp = 0) is then 


W 

LPTLS 


2 


1 (nd)p(o) 


(55) 


Figure 7 clearly reveals that the initial state has to be 
carefully prepared in order to exhibit the maximal possi¬ 
ble power boost. 


VII. DICKE-SYSTEM PERFORMANCE 


with an arbitrary phase $. Alternatively, these condi¬ 
tions hold (since the coherences between the ground- and 
excited states do not matter) for the incoherent mixture 

P(0) = Poo(O) |0)(0| + [1 - poo(O)] (53) 

in the new basis containing the bright state ]cf. Eq. (26)] 

= ( 54 ) 

V1 + 

Any deviation from this “optimal” initial condition results 
in a reduced thermalization capability and therefore in a 


Let us now consider an ensemble of N two-level atoms 
interacting with the same bath. The bath mediates the 
dipole-dipole interaction as well as non-local cooperative 
photon exchange, whereby a photon emitted by one atom 
can be re-absorbed by another ]42, 43]. Quite gener¬ 
ally, the dipole-dipole interaction renders the j-excitation 
states non-degenerate ]43-45] . Whilst some of these col¬ 
lective states give rise to an enhanced decay rate, others 
decay much slower than their single-atom counterparts 
(see Fig. 8), conforming to “superradiance” and “subra¬ 
diance”, respectively ]42-52]. Hence, in general V-atom 
systems do not obey the Dicke model ]22]. However, 
by choosing an appropriate geometry (elaborated on in 
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Figure 8. (Color online) The Dicke system exemplified for 
N — 2. The two atoms interacting with the same environ¬ 
ment (left plot) are mapped onto a non-degenerate four-level 
system consisting of collective states (middle plot). The two 
states carrying one excitation are superradiant (|1)) and sub¬ 
radiant (I l^)) respectively. If the dipole-dipole interaction be¬ 
tween the atoms vanishes (fli 2 = 0) and their collective decay 
rates are maximized (ri 2 = 7 ), the system is mapped onto an 
effective three-level system with equidistant levels (equivalent 
to a spin-1 system plus one dark state) (right plot). 


Sec. Vine), the dipole-dipole interaction can be made 
to vanish, while at the same time the subradiant states 
become completely dark (see Fig. 8 ). Only such a config¬ 
uration can realize the Dicke system —N two-level atoms 
(equivalent to spin- 1/2 particles) can then be mapped 
onto a collective spin-A^/2 system [25], whose Hilbert 
space is spanned by lV-l-1 (symmetrized) collective states 
\j) containing j excitations with eigenenergy jhujQ, where 
ujQ denotes the TLS transition frequency, and 2^ —(iV-|-l) 
dark states, which are decoupled from the dynamics (By 
contrast, in the fully aligned multilevel case N — 2 states 
of the total N states are dark). This mapping [25, 43[ is 
sketched in Fig. 8 for two atoms. 

The dissipative dynamics of the Dicke system is de¬ 
scribed by the Liouvillian superoperator (for brevity 
given here for a single bath and without modulation) 
(see also Ref. [25[) 


1 


Cp = N-G{ujo) {2ApA"^ — A^Ap — 

+ {2A^pA - AA^p - pAA^) 


(56) 


with 


1 - (Hd) 


p(o) 


agonal, partially-thermalized, state 
L 

J = l,...,iV 


(58a) 

(58b) 


where we again resort to the projector Hd onto the dark 
subspace. Here the states |(/j) again denote the eigen¬ 
states of the dark part of the initial density matrix. The 
major differences of Eqs. (58) to Eqs. (36) are the Boltz¬ 
mann factors, which now differ according to the number 
of excitations in each state |j). As in Eqs. (36), coher¬ 
ences between the bare atomic states explicitly appear 
when transforming the solution (58) to the 7V-atom prod¬ 
uct basis. 

The same derivation as before for the multilevel case 
yields the cold and hot heat currents presented in Ap¬ 
pendix F. The ratio of these heat currents (FI) and 
the associated power to their counterparts generated by 
N independent two-level heat machines [Eq. (El) with 
N = Aeff = 2[ reads 


J^ 


W 




N-l 

g-jPelfULJo 

j=o 


^ (59) 

tls ’ 

Poo 


where 


Poo 






1 - (Hd) 


p(o) 


(60) 


At low effective temperatures, /3eff —>■ 00 , at most a 
single excitation is present in the Dicke system. In this 
limit, the Dicke ladder for N — 1 atoms may be mapped 
onto a degenerate TV-level system with parallel transition- 
dipole vectors (which corresponds to an effective two- 
level system with transition-dipole strength enhanced by 
y/N — 1). Now, from Fig. 6 we know that in the low- 
temperature regime the aligned and non-aligned multi¬ 
level systems perform similarly, the latter corresponding 
to — 1 independent TLS. This behavior is obtainable 
also from Eq. (59), 


W 


lim 


1 - (Hd) 


p(o) 


(61) 


N N-l 

i=l j=0 

and the individual Pauli operators := \gi){ei\ for the 
ith two-level atom. The jump operator oc y/NjA de¬ 
scribes the emission of one quantum accompanied by de¬ 
excitations (j —>■ / — 1 ) of the non-dark excited states 
\j > 0). The individual emission rates from \j) to \j — 1) 
(j S (0, A^]) evaluate to = 'yj{N — / + 1) [25[. The 
steady-state solution of the master equation with the Li¬ 
ouvillian (56) is then—in analogy to Eqs. (36)— the di- 


Hence, at low effective temperatures the N-atom Dicke 
heat machine performs as well as (but not better than) 
N independent two-level atoms, provided it is amenable 
to complete thermalization. 

In the opposite limit of high effective temperatures, 
PeS —>■ 0, the cooperative Dicke system yields a power 
boost. 


W 


lim 

/3eff-*-0 A^VFTLS 


2N 
N + 1 


1 (nd)p(o) 


(62) 


This ratio has the following meaning: The denominator 
represents the thermal equipartition of all non-dark state 






































13 


.9 

’■+J 

^-1 


O 


% 

o 

a, 



Figure 9. (Color online) Heat currents and power output 
[Eq. (59)] for the A'^-atom Dicke model compared to N in¬ 
dependent two-level-atom heat machines for different particle 
numbers. Here we have assumed an optimal initial condition 
orthogonal to the dark subspace. The maximum enhancement 
factor is 2. 


| 2 . 0 , 0 ) | 2 , 1 ,± 1 ) | 2 , 1 , 0 > 


12 , 0,01 + 12 . 1 . 0 ) | 2 , 1 .± 1 ) 
■/i 


| 2 , 0 , 0 )-| 2 , 1 , 0 ) 

^/2 


electric field 


11)0,0) 


|1)0,0) 


Figure 10. (Color online) Linear Stark effect in hydrogen. 
The metastable state |2,0,0) is mixed with |2,1,0) to form 
the dressed states (|2,0,0) ± |2, l,0))/-\/2 by a static electric 
field [53]; for a weak field these states are quasi-degenerate. 
Their transition dipoles (red arrows) are anti-parallel, but 
since |2,0,0) and |2,1,0) have opposite parity, a hydrogen 
atom in a static electric field behaves as if it possesses a per¬ 
manent dipole moment ].53]. 


VIII. REALIZATION CONSIDERATIONS 
A. Misaligned dipoles 



maximal power 
relative to TLS 

(low Tefl) 

maximal power 
relative to TLS 
(high Tefi) 

initial state 
for maximal 

power 

A-level 

(non-aligned 

dipoles) 

N -1 


any 

A-level 
(fully aligned 
dipoles) 

N-1 

N-1 

non-dark 

A-atom 
Dicke 

N 

2N^ 

N+1 

non-dark 


Table 1. Summary of the (maximal) power enhancement fac¬ 
tors and their respective dependencies on the effective tem¬ 
perature ]cf. Figs. 6 and 9] and the initial condition for the 
multilevel heat machine and the Dicke machine (compared to 
a single TLS). The Dicke heat machine shows the best per¬ 
formance. 


populations. 


1 — (nd)^(Q) ]/(IV -I- 1), whereas N is the 


bright-state power enhancement and the factor of 2 stems 
from the high-temperatnre TLS ground-state population 
1 / 2 . 

Hence, the iV-atom Dicke heat engine can give at most 
twice the power of its connterpart consisting of N in¬ 
dependent TLS. The fnll temperature dependence of the 
power boost for varions ensemble sizes is shown in Fig. 9. 

As a consequence of Eq. (59), i.e., the identical boost of 
the two heat currents, the efficiency of the Dicke heat en¬ 
gine and the COP of the Dicke refrigerator, respectively, 
are the same as their two-level counterparts. 


^Dicke ^ ^TLS 

COpDick® = COP^^® 


(63a) 

(63b) 


Onr results on power enhancement in iV-level systems 
as well as for N atoms conforming to the Dicke model 
are snmmarized in Table 1. 


Degenerate excited states with non-parallel dipole- 
transition vectors to the ground state are nbiquitous in 
atomic systems withont hyperfine splitting. For exam¬ 
ple, the 2p manifold of hydrogen consists of the three 
degenerate states |n = 2,Z = 1 ,to = {0,±1}). Their cor¬ 
responding dipole transition vectors to the gronnd state 
are all of eqnal strength, but orthogonal to each other. 
This orthogonality, however, does not restrict thermal- 
ization (bnt may not produce the maximal boost at high 
temperatures, cf. Fig. 6), since there is no difference be¬ 
tween orthogonal and any other non-aligned confignra- 
tion in the steady-state heat machine. 


B. Aligned dipoles in multilevel atoms 

Multilevel degeneracy combined with perfect align¬ 
ment (parallelism) of the transition dipoles, which is a 
prerequisite for coherence effects, is much harder to re¬ 
alize, since it cannot occur in bare atomic systems due 
to selection rules [38]. It is, in general, also not possible 
to engineer dressed states within the excited manifold 
that have this property and at the same time adhere to 
our model: Combining degenerate states with orthogo¬ 
nal dipole moments (as in the 2p manifold of hydrogen) 
creates effective states with non-orthogonal but also non¬ 
parallel dipoles [38]. 

A possible realization of alignment may be to super¬ 
pose a decaying and a non-decaying (metastable) state. 
Such linear combinations of decaying and metastable 
states, however, do not possess a definite parity. This 
means that the dipole operator D has permanent (static) 
dipole moments (ei|D|ei), as e.g., formed via the linear 
Stark effect in hydrogen (see Fig. 10). In the derivation 
of the master equation (15a), however, we have assumed 
that the ground and excited states have a definite parity 
and therefore the atom has no permanent dipole moment. 

Nevertheless, a permanent dipole moment (ei|D|ei) is 
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not a hindrance to the treatment of system-bath cou¬ 
pling for a broad class of thermal baths [25]: For dipolar 
system-bath coupling [Hsb = D • B as in Eq. (1)], the 
static part of the dipole operator in the interaction pic¬ 
ture does not contribute to the Liouvillian part of the 
master equation as long as 

lim G{ui) = lim "f(u}){n{uj) -I- 1) = 0, (64) 

cj—vO cj—>-0 

which is the case for bosonic baths provided 7(0;) oc a;“ 
with a > 1. This uj = 0 (static) component of the dipole 
operator, which is absent in the Liouvillian, represents 
the permanent dipole moment, since (given here for a 
TLS for brevity) 

D(t) = da+e“«‘ -h -f (e|D|e) a+a. (65) 

with the transition dipole matrix element d := (e|D|g). 

Yet, even though Stark-shifted states may realize the 
aligned-dipoles model (on time scales much shorter than 
their inverse splitting 1 /A), the predicted power boost is 
countered by the fact that in an equal superposition of de¬ 
caying and metastable states only half of the superposed 
state has dipolar interaction with the bath [21]. Hence, 
the decay rate is reduced compared to that of the decay¬ 
ing state to 7 eff = I 7 . Consequently, the power boost 
and the reduction in the effective spontaneous emission 
rate exactly cancel each other. Still, one could argue, 
that there is a power boost relative to a TLS decaying at 
rate 7 eff. 

It has been suggested that one may create an effective 
E-type system with aligned dipoles out of a A-system 
with orthogonal transitions by applying a strong laser 
field [38] . The resulting E-system is non-degenerate, but 
the energy mismatch can be tuned to be very small, such 
that one can observe the effects discussed here on long 
time scales. However, as for the Stark effect discussed 
above, the decay rate of the dressed state is smaller than 
that of the bare excited state, precluding a power boost. 

In view of these difficulties, the strict alignment re¬ 
quirement should be relaxed in realistic experimen¬ 
tal situations. Even though for quasi-aligned dipole 
transition vectors no dark states exist, strongly sub¬ 
radiant (metastable) states may still arise. The non- 
thermalization of these “quasi-dark” states should there¬ 
fore be experimentally observable on appropriate time 
scales. 

A promising alternative is to make use of the plethora 
of vibrational and rotational levels in molecules as dis¬ 
cussed in the supplemental material of Ref. [11]. 

C. The Dicke system 

In view of the conceptual difficulties in finding a real¬ 
ization of parallel dipoles in a degenerate multilevel sys¬ 
tem that is capable of power boost compared to a TLS, 
we shall revisit the Dicke system discussed in Sec. VH. 
As we have seen, this multipartite setup exhibits similar 




Figure 11. (Color online) A realization of the Dicke model by 
means of two-level atoms confined at interatomic distance A 
next to or within one-dimensional waveguides. See text for 
details. 

thermodynamic properties to a degenerate Y-level sys¬ 
tem of parallel dipoles and may produce a power boost 
under suitable initial conditions. 

Traditionally, the superradiant Dicke model was 
thought to be realizable by a dense atomic ensemble con¬ 
fined well within the cubed emission wavelength [25, 35] . 
Yet, in general, such ensembles do not obey the Dicke 
model, since the dipole-dipole interaction as well as the 
cooperative decay rates strongly vary, depending on the 
spatial symmetry and spacings of the atomic ensem¬ 
ble [42, 54-58]. The remedy is to realize the Dicke model 
in field-confining structures such as photonic bandgap 
structures [59, 60], cavities [61] (see also [62]), and waveg¬ 
uides [63-66]. 

In particular, consider the realistic setup (Fig. 11) [67] 
of a periodic Id lattice of atoms (each tightly confined 
by the lattice potential) that are positioned at equal dis¬ 
tances d within or next to a photonic waveguide. Al¬ 
ternatively, we may employ a chain of superconducting 
qubits coupled to a microwave coplanar waveguide [ 68 ] . 
Field confinement can give rise to a giantly enhanced 
density of modes at the lower cutoff frequency of the 
waveguide [64, 69] . Then, atoms whose resonance fre¬ 
quency ojQ is just above Wcutoff are predominantly cou¬ 
pled to the Id axial-mode continuum, whereas transverse¬ 
mode (free-space) effects are negligible in comparison, 
7id 3> 7free [64, 69]. In such an effective Id photonic 
“bath”, the cooperative decay rates Fy and resonant 
dipole-dipole interactions of atoms i and j are [ 68 , 70] 

Tij = 7id cos (kd) (66a) 

ilij = sin (kd ), ( 66 b) 

with the wavenumber k = 27: j\ = loq jc. The oscillatory 
behavior of Eqs. ( 66 ) stands in stark contrast to the decay 
with distance of their counterparts in free space [42] and 
allows us to choose the interatomic distance d = X such 
that 

Tij = 7 id (67a) 

= 0. (67b) 

In this geometry, the atomic lattice realizes the Dicke sys¬ 
tem [22, 35, 44], i.e., diagonalization of the Y-atom Li- 
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ouvillian operator with the interatomic decay rates (67a) 
leads to the Liouvillian (56) involving a single decay chan¬ 
nel, whose collective steady-state solution is the partial 
thermal state (58). By coupling the atoms to two similar 
waveguides kept at different temperatures (see Fig. 11), 
we may operate the setup as the heat machine depicted 
in Fig. la. 

The high-temperature regime, where the Dicke model 
may outperform the fully thermalized model corresponds 
to temperatures Teff > huJo/{5kB) (see Fig. 9). For opti¬ 
cal transitions (wq 10^^ Hz) these are temperatures of 
the order of 10^ K. By contrast, for microwave transitions 
(wn 10® Hz) the high-temperature regime is already re¬ 
alized at Teff > 10-3 K. 


IX. CONCLUSIONS 

The present analysis sheds new light on multilevel or 
multipartite heat machines and their degeneracy as a 
thermodynamic resource that can boost the heat cur¬ 
rents and the power output of heat engines. By contrast, 
the efficiency of such machines was shown to be unaltered 
compared to a TLS-based machine, which adheres to the 
Carnot bound. 

We can summarize our findings regarding the role of 
initial and steady-state coherences in the power enhance¬ 
ment of degenerate steady-state multilevel quantum heat 
engines as follows, 

1. Steady-state coherences are a consequence of ther- 
malization if at least two transition dipoles are par¬ 
allel. Such coherences arise in the bare energy basis 
but not in the rotated basis containing dark states. 

2. Neither initial nor steady-state coherences do neces¬ 
sarily imply power enhancement compared to a sin¬ 
gle TLS and, vice versa, power enhancement does 
not always entail initial or steady-state coherences. 
The relevant factor for power enhancement is the 
thermalization capability (Sec. VB) of the initial 
state, not the amount of steady-state coherence. 
As an extreme example, even if all three dipoles of 
a four-level system are perpendicular to each other 
(the largest orthogonal system possible in three- 
dimensional space) the power is enhanced—without 
coherence between the excited states at any time 
and regardless of the initial condition. The decisive 
ingredient is the common ground state of the dipole 
transitions. It is generally true that aligned dipole 
configurations yield higher power boost than non- 
aligned ones in the high-temperature limit—if the 
initial state is properly chosen. Even then, we may 
not attribute this power boost (compared to non- 
aligned configurations) to steady-state coherences 
but rather to the presence of dark states, which 
reduce the number of available levels amenable to 
thermalization. The key factor for enhancement is 


the ground-state population of this reduced Wff- 
level system. The fewer levels can thermalize, the 
more is the ground state thermally populated. Yet, 
any population in the ground state necessarily de¬ 
creases the maximally permitted coherences within 
the excited-states manifold. 

3. Power reduction implies the existence of steady- 
state coherences. Power reduction is associated 
with dark states being initially populated. The 
largest possible steady-state coherences correspond 
to fully excited (population-inverted) states, which 
can only occur if the initial state is dark. 

Overall, we found that an Wlevel heat machine can¬ 
not yield higher power than TV — 1 independent two-level 
heat machines. Whilst at low-temperatures both ma¬ 
chines perform equally well, the power enhancement in 
the iV-level machine decreases at higher temperatures, 
unless all dipoles are parallel. Taking advantage of the 
full cooperativity of the excited states (occurring only for 
parallel transition dipoles), however, requires carefully 
chosen (non-dark) initial conditions (Sec. V) and its re¬ 
alization encounters conceptual difficulties (Sec. VHIB). 

Similar conclusions were found to apply in Watom 
ensembles that, under specially-chosen conditions, obey 
the Dicke model (Sec. VHIC). Remarkably, the Y-atom 
Dicke model yields at best doubly-enhanced power com¬ 
pared to N independent atoms (and only at high tem¬ 
peratures) (Sec. VH). Nevertheless, doubling the power 
output of a heat machine by cooperative effects will be 
an important achievement. 

To conclude, our findings elucidate the possible im¬ 
pact of steady-state coherence or entanglement in the 
working medium on the performance of heat machines. 
The fact that the heat baths are assumed to be thermal 
allows us to uphold the Carnot bound. The possibility 
of exceeding the Carnot bound by virtue of harnessing 
non-thermal baths in the heat machine, including baths 
consisting of coherently-superposed [7, 71, 72] or superra¬ 
diant [73] atomic systems, as well as squeezed baths [74] 
is outside the scope of the present discussion, and so are 
effects related to piston quantization [75-79] or to strong 
system-bath coupling [80-82] that require further clari¬ 
fication. Notwithstanding these open issues, the present 
results help delineate the limited (but significant) part 
that quantum coherence or entanglement may play in 
quantum thermodynamics. 
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Appendix A: Coefficients and solution of the ODE 
for the three-level system 


Appendix B: Proof of the steady-state solution (24) 
for non-aligned dipoles 


The coefficient matrix A and the inhomogeneity b ap- Proof. The Liouvillian (15b) for a single heat bath at 


pearing in the ODE (18) are (generalizing Ref. [21]) 


A -G(a;o)x 


inverse temperature /3 and without modulation (q = 0) 
reads 


Cp = Ig{uo) E 


N-1 


/ —1 — a 


2 0 pae^'f’{1 + 26-^^^^°) 0 

0 -1-a^ pae-^^ {l + 2e-^^‘^°) 0 


\ 


2pae 


— lip 


2pae*‘^ -2 - (1 + a^) 2{a^ - 1) 


y—pae 


— hnip — 


pae 


Ip 


and 


b := -G(a;o) 


2Q,2g-/3fii^o 


-pae*'^ 
—pae”®"^ 


-2a^ / 


(Al) 


i=i 

N-1 


+ E 




+ -G(-a;o) 5] 


a 


PV + i: fn 


>e ^33'afa 




oAo 


with the abbreviations T)]^ := V 


(A2) 


V 




i 




(Bl) 


and := 


L ). We now apply this operator to the diagonal 


ansatz 


N 


Here we have defined the shorthand notation pi 2 = 
p 2 i =: p, pi 2 = -<<321 =: ^ [cf. Eq. (4)], and the dipole 
strengths are ai = 1 and a 2 ='■ a [cf. Eq. (3)[. 

The coefficient matrix appearing in the ODE (20) ex¬ 
plicitly reads 




(B2) 


k=l 


a 0 c* ON 
0 a 0 c 


where 


c 0 I 0 

vO c* 0 6 


a := 1 + (1 + af 

5 := + (1 + a)2 e-/3R-o 

c := pae®'^. 


(A3) 


(A4a) 

(A4b) 

(A4c) 


Making use of the KMS relation (16) we find 
= iG(cuo)x 

N-1 

X E {‘^PU - (|0)(0| - K)(z|) - 


2=1 




N-1 


X E (-P- - p- + 2 plle-^^-«) \i){j\. 

(B3) 


ij'=l 


Hence the diagonal thermal state pf^ = e 
The general steady-state solution of the linear [Eq. (24)[ is indeed the steady-state solution of the mas- 


ODE (18) is then found to be 

.iv l-l(2-|-e'5'“"0) 


II 

' 1+a^ 

[0 

pr2 = (p?i)* 

1 

f 1 (IT') 

Poo = \ 

1 i+e-Ohu^o fJ- -k) 

1 1 

1 

1, l+2e-^^‘^o 
^ 1 a^+X(l — 

II 

1 1+q 2 i_|_g/3»ii^0 

1 1 


P = 1 
P e [0,1) 

p = i 

P e [0,1) 


2+e/3'*‘"o 


P = 1 

P e [0,1) 


(A5a) 

(A5b) 

(A5c) 

(A5d) 


ter equation £p®® = 0 for A levels when the dipoles are 
not aligned. 


Appendix C: Proof of the steady-state solution (28) 
for aligned dipoles 


Proof. The master equation 

P — (Gdecay “f ^absorption) P 


(Cl) 


associated to the system-bath interaction Hamilto¬ 
nian (25) is the sum of the two Liouvillians 


where 


'IV-l 


1 


GdecayP= xG(wo) E “j ^(o’-W+) (C2a) 


21{®^Pii(0) + P22(0) — 2a Re [e®‘^pi2(0)] } 

(A6) 

is an integral of motion associated with the overlap of the 
initial state with a dark state, to be discussed in Sec. V. 


. Al 

^ N-1 


GabsorptionP = -G(-a;o) ( E “i ) (C2b) 

. Al 
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describing decay and absorption processes mediated by 
the bath between the bright- and the ground states, re¬ 
spectively. 

We first show that the steady-state solution cannot 
contain coherences between non-dark and dark states be¬ 
cause according to equation (C2) the Ehrenfest equations 
for these coherences read 

(|V'd')(V’b|) oc-(|?/'d)(V'b|) (C3a) 

^(|^,^)(0|)cx-e-^"-“(|V^d")(0|), (C3b) 

and hence (|V'd)(V^b|),,,, = 0 and (|V’d)(0|)ss = O' 
sequently, the steady-state solution can be decomposed 
into a dark-state portion and a portion in the subspace 
orthogonal to the dark states, 

_ss ss I _ss f A\ 

P = P±+ Pd (C4) 

with the normalizations 

TrpJ = 1 - (nd)^(o) = (V'blp(O)IV'b) + (0|p(0)|0) (C5a) 

Trp- = (nd)^(o) ^ ^ 
i=i 

The latter normalization condition follows from the fact 
that the relative weights of non-dark (bright and ground) 
states and initial dark-state contributions cannot be al¬ 
tered by the Liouvillian dynamics. 

We then assume a diagonal steady-state density ma¬ 
trix. Inserting such an ansatz into the master equa¬ 
tion (C2) and making use of the KMS relation (16) yields, 
in the basis spanned by the bright and ground states. 


strength between the ground state and a certain bright 
state. Accordingly, the interaction of the Zth “domain” 
(1 < ^ < p) with the bath is described by the effective 
Hamiltonian 


jgtI _ 

■^SB, domain a 

ni 


\ 

j=ni-i + l 

X (d(|_ (g) e/ • B + ctL (8) e* • B^) 

(Dl) 

with the Pauli operators [cf. Eq. (27)] 



cr+ := V'b)(0 

(D2a) 


ctL := 0)(V'b 

(D2b) 


and the unit vector e; G pointing along the transition- 
dipole direction of the Ith domain. 

The interaction of the remaining transitions, which are 
not parallel to any other, with the bath are described by 
the Hamiltonians {Np -1- 1 < Z < — 1) 

H^Sb = azldil {a^^ei-B + ai^e* ■ B^) , (D3) 

with the “bare” Pauli operators := |0(0| and aL ■= 

|0)(/| defined in Eq. (2), and e; pointing in the direction 
of the ffh dipole transition vector. 

The interaction Hamiltonian can thus be decomposed 
as 

p N—l 

HsB = H^B, domain + ^ ^gB' (^4) 

;=i i=Np+i 

It formally corresponds to a system with 


P(0) V'd) • 


(C5b) 




|0)(0|-|^b)(^b| . (C6) 


The right-hand side obviously gives zero for = 

p^Qe~^^°. Combining the latter relation and the nor¬ 
malization condition (C5a) yields the populations p^^ 
and Pqq from equation (28). 

The part of the initial density matrix lying within the 
dark-state space is not altered by the dynamics and hence 
= Pd(0). This state, in principle, can be an arbitrary 
superposition of the dark basis vectors. We aim at a 
diagonal representation of the full density matrix, hence 
we diagonalize the dark-state density matrix and find the 
relative weights which are the corresponding eigen¬ 
values of pd(0). ■ 


Appendix D: System-bath interaction Hamiltonian 
for arbitrary dipole geometry 

As shown in Sec. VB, parallel dipoles effectively corre¬ 
spond to a single transition with an enhanced transition 


N,g:=p + N- Np (D5) 

levels that consist of the p bright states, the — 1 — Np 
excited states with non-aligned transition dipoles, and 
the ground state |0). The remaining ~ 1) = 

N — Nes levels are dark and thus decoupled from the 
dynamics. 


Appendix E: Heat currents and power for the 
multilevel system 


The heat currents associated to the Ah bath (i € 
{c, h}) for the degenerate A-level system explicitly read 


Ji = ^ HuJo + qVt)P{q)Gi{ijjQ + qil) X 
q& 

^-Pih(LJo+qQ,) _ I ^ 

1 + (A^eff ~ 1 ^ 

(El) 

Having found the heat currents, we may now compute 
the power, which according to the first law of thermo¬ 
dynamics (energy conservation) reads W = —{Jh + Jc)- 


j 1 (nd)p(o) 
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Upon inserting the heat currents (El) for i G {c,h} into 
this equation, we find the power 


W = E E h{uJo + qft)P{q)Gi{uJo + gU)x 

qGZ iG{c,h} 

_ ^-/3ih{uJo+qQ) 

1 + (A^eff ~ l)e“'®»ft^‘^o 

(E2) 


' N-l 

E« 


1 - (Hd) 


p(o) 


We again observe the dependence of the heat flows on 
the thermalization capability 1 — (nd)p(Q). These heat 
currents are governed by three inverse temperatures, /3c, 
^h, and ^eff- 


Appendix F: Cooperative heat currents for the 
Dicke system 

For the 7V-atom Dicke system the heat currents read 
{i G {c,h}) 


Ji — 'y ' Nti(ujQ+qil)P{q)Gi{uJo+q^l) 

q^Z 


N-l 


3=0 


X 


X 


^-l3ih(uo+qCl) _ g-/3etfS.<^0 

p-iPe!ihuio 
l^j=0 ^ 


1 (nd)p(o) 
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